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The long-time asymptotic solution of the Korteweg-de Vries equation for general, step-like initial data is analyzed. Each sub-step in 
well-separated, multi-step data forms its own single dispersive shock wave (DSW); at intermediate times these DSWs interact and 
develop multiphase dynamics. Using the inverse scattering transform and matched-asymptotic analysis it is shown that the DSWs 
merge to form a single-phase DSW, which is the 'largest' one possible for the boundary data. This is similar to interacting viscous 
shock waves (VSW) that are modeled with Burgers' equation, where only the single, largest-possible VSW remains after a long 
time. 
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1. Introduction 

Dispersive shock waves (DSWs) appear when dispersion 
dominates dissipation for step-like data; they have been seen 
in plasmas [1], fluids (e.g., undular bores) El [3], superfiuids 
El E 0, and optics E [TUl EJ. The Korteweg-de Vries 
(KdV) equation is the leading-order asymptotic equation for 
weakly dispersive and weakly nonlinear systems Ifl2l . Each 
step in well-separated, multi-step data forms its own DSW; 
these DSWs interact and develop multiphase dynamics [13|. 
Here we show that these DSWs merge in the long-time limit to 
form a single-phase DSW; the boundary data determine its form 
and the initial data determine its location. We find this long- 
time asymptotic solution using the inverse scattering transform 
(1ST) and matched-asymptotic analysis. This merging of shock 
waves is similar to viscous shock waves (VSWs), where dissi- 
pation dominates dispersion and only a single VSW remains in 
the long-time limit. While the linear KdV equation's solution 
for step-like data in the middle region has width 0{t 1 ^) (see 
{j |4.2| i and the KdV equation's solution for vanishing data has 
a collisionless shock of width (9[f^ 3 (logf) 2 ^ 3 ] (see 5 4.3 1, the 
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Figure 1 : Numerical solutions of the KdV equation for the step and the van- 
ishing initial data shown in gray, (a) For step data, there are three basic re- 
gions: exponential decay in region A, the DSW in region B with width O(t) and 
height 0(1), and an oscillating tail in region C. (b) For vanishing data, there 
are four basic regions (see 1181 ). The collisonless shock in region III, which 
is analogous to the DSW in region B, has width <9[r l/3 (log r) 2/3 ] and height 
<9[(log0 1/2 r 2/3 ]. 



u(x, f) goes rapidly to the boundary conditions 



lim u — and lim u = -6c , 

X — > — OO A" — > + 00 



(2) 



DSW that we find has width 0(t) ($32). We anticipate that our 
1ST and matched-asymptotic procedure will be applied to other 
integrable nonlinear partial differential equations (PDEs) with 
step-like data. 

Here we consider DSWs described by the KdV equation — 
the leading-order asymptotic equation for systems with weak 
dispersion and weak, quadratic nonlinearity. The KdV equa- 
tion, in dimensionless form, is 



u t + uu x + s u x 



0. 



(1) 



where subscripts denote partial derivatives. We require that u 
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where s and c are real, positive constants. Here, e corresponds 
to the size of the regularizing dispersive effects. (Since the KdV 
equation is Galilean invariant, we can transform any boundary 
conditions where lim^-co u > lim A ^ +00 u to these boundary 
conditions.) We use the 1ST method (see OH [15] [RO [HI) and 
matched-asymptotic expansions (see lfl8l [191 ) to find a large- 
time asymptotic solution. 

Single-step data, such as a Heaviside function, evolve to form 
a single DSW. This DSW has three basic regions (Fig. [lj): a 
rapidly decaying region to the right of the DSW (region A); 
the central DSW (region B), which is a slowly varying cnoidal 
wave with a soliton train on its right and an oscillatory tail on its 
left; and a decaying, oscillatory region to the left of the DSW 
(region C). Similarly, each step in well-separated, multi-step 
data forms its own DSW (Fig. |2^); these DSWs eventually in- 
teract strongly and at intermediate times exhibit multiphase dy- 
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Figure 2: A numerical simulation (using the scheme in 1131 ) of three well- 
separated steps (at ; = 0) that: (a) form three single-phase DSWs at f = 1 ; (b) 
and (c) interact strongly and exhibit multiphase dynamics; and (d) eventually 
merge to form a single-phase DSW. Here, e 2 =0.1 and c = 1. 



Hilbert problems that depend on a parameter (here time, t ); they 
have been used to find the long-time asymptotic solution for 
vanishing data (see BT1 [32)) — see also l33l for a nonlinear 
Schrodinger equation shock example. For our purposes here, 
the matched-asymptotic method is sufficient. 

In this letter: We give the required 1ST results in ^2] Then 
we find the rapidly decaying solution in the region to the right 
of the DSW (pi). This matches into the central DSW (p^2), 
which is a slowly varying cnoidal wave with a soliton train on 
its right and an oscillatory tail on its left. Finally, we match the 



decaying, oscillatory region to the left {{ 3.3 i with the DSW. We 
then compare this solution with: the solution of Burgers' equa- 



tion (54.1 1, which is the leading-order asymptotic equation for 
VSWs; the solution of the linear KdV equation with step-like 



data (§4.2 1; and the solution of the KdV equation with vanish- 



ing data (§4.3 1. Then we draw some conclusions (54.4i 



2. 1ST solution 



namics (Figs, and|2j;). In this letter we show, in the long- 
time limit, that these DSWs eventually merge to form a single- 
phase DSW (Fig. [2}I). 

We use 1ST and matched-asymptotic methods to show that 
general, step-like data go to a single-phase DSW in the large- 
time limit. Individual, single-phase DSWs have been exten- 
sively studied using wave averaging techniques (see [20, 21 , 22 
23 1), often referred to as Whitham theory [24, 25 1. The evolu- 
tion of two-phase DSWs to a single-phase DSW was investi- 
gated by Grava and Tian ll26l in the zero-dispersion (e — » 0) 
finite-time limit using Whitham theory and by Ablowitz et al. 
|[T3l in the fixed-dispersion long-time limit using numerical and 
asymptotic methods. Both zero-dispersion and large-time are 
important, but different, limits; here we study the large-time 
limit with fixed dispersion. By using the 1ST method we find 
the asymptotic solution directly; in Whitham theory, we must 
evolve the solution through intermediate times. Therefore, we 
can investigate general, step-like initial data and the interaction 
of DSWs without having to find the solution at intermediate 
times. 

The 1ST theory for step-like initial conditions was studied 
in E3 and J28] |29l. Hruslov E3, based on (30l, gives the 
Gel'fand-Levitan-Marchenko (GLM) integral equations and 
investigates the soliton train associated with the DSW. Co- 
hen [28] and Cohen and Kappeler |29l , using the methods in 
lfT31[30l . rigorously studied the properties of the scattering data, 
rederived the GLM integral equations, and analyzed existence 
of solutions corresponding to certain initial conditions. We state 
the 1ST results that we need to find our asymptotic solution in 

From these 1ST results, we use and suitably modify the 
methods in lM~8l [T9l to find our long-time asymptotic solution. 
Ablowitz and Segur lfT8l [T9l developed the 1ST and matched- 
asymptotic methods for vanishing data (where c = 0); we mod- 
ify them for step-like data (where c + 0). The results of the 
long-time asymptotic analysis in this letter are new. There are 
elegant and powerful asymptotic methods based on Riemann- 



The 1ST method transforms the initial data into scattering 
data, evolves the scattering data in time, and then recovers the 
solution from the evolved scattering data. First, we associate a 
linear (Lax) pair with the nonlinear PDE, in this case ([1]). Then 
we use the scattering equation of the linear pair to transform 
the initial data into scattering data. The scattering data are then 
evolved in time using the associated linear equation. Finally, 
the solution is recovered using a linear integral equation, the 
GLM integral equation, at any time. 

The Lax pair associated with ([TJ is 

V xx + (u/6 + A 2 ) v/s 2 = 0, (3a) 
v, = (u x /6 +y)v + (4A 2 - u/3) v x , (3b) 

where y is a constant. This linear pair is compatible (that is, 
v xx t = Vtxx) when u satisfies ([T| and A is isospectral (that is, 
dA/dt = 0). The eigenfunctions that satisfy ( f3a| > are defined 
using Q: 

4>(x; A) ~ e~ iAxlB , 4>(x; A) ~ e iAx/B , (4a) 

as x — > -oo and 

tfr(x; A r ) ~ e a ' x/s , if/(x; A r ) ~ e ~ a ' xls , (4b) 

as x — > +oo, where A r = V/l 2 - c 2 . The branch cut of A r is taken 
to be A e [-c, c]; the branch cut of A is taken to be A r e [~z'c, ic\. 
so Im(/l r ) ^ when \m(A) ^ 0. This branch cut is one of the 
main differences between vanishing and step-like data. 

The Wronskians W(4>,tf>) = liA/e and W(i//,if>) = -2iA r /e 
are constant. The scattering eigenfunctions and scattering data 
a and b associated with ( f3a| > satisfy 

4>(x; A) = a(A, A r )\j/{x; A r ) + b(A, A r )if/(x; A r ) (5) 

for A r + 0, A r e K (or, equivalently, \A\ > c, A e R). The 
scattering data can be written as 2iA r a = eW(<p, if/) and 2iA r b = 
sW(ip,<p). We can use this to extend a to A e (-c, c), where A r 



2 



is pure imaginary; when A e (— c, c), if/ is real and exponentially 
decaying. 

We define the transmission coefficient T = l/a and the re- 
flection coefficient R = b/a so that |5]) can be written as 

T(A, A r )<f>(x; A) = \fr{x; A r ) + R(A, A r )i//(x; A r ). (6) 

From (|2j, pb] i, and (|6j, these coefficients evolve in time as 
T(A, A r ; t) = ^ 7XA, A r ; Qy^^+^Ms and ^ ^. f) _ 
fl(/l, /l r ; 0)e /(8/l2/l ' + 4c2 'W/ £ . Unlike with vanishing data, the trans- 
mission coefficient, T, depends on time; this time dependence 
is not purely phase when A e (-c, c). 
The associated GLM integral equation is 

Xoo 
Q.iy + z;t)K{x, Z ;t)dz = 0, (7) 

where 



i r°° 

^ ;t) = 2^h ReUrS/EUr + ^ Cje 
+ 2sn X 



WIArte 



dA, 



the constants \iKj} ! J =l are the (simple) poles of T(iKj,A r (iKj); t), 

Kj = tJk 2 : + c 2 , cj — -i/ij/(sdA r a(iKj)), (f>(x; iKj, t) = 
fj.j(t)if/(x; iKj, f), and < k\ < ••• < kh are real. We omit any 
contributions from poles in our asymptotics: the poles are re- 
lated to the solitons, which move to the right, and so do not 
affect the DSW in the large-time limit. 
From K, we recover u(x, t) from 



u{x,t) = -6c 2 + 12e 2 — K(x,x;t). 

dx 



3. Long-time asymptotics 

To find the long-time asymptotic solution: we find the so- 
lution right of the DSW, then we use it to find the DSW, and 
then we find the solution left of the DSW and match it into 
the DSW. We use |7]i to asymptotically compute the behavior 
right of the DSW. When this asymptotic solution breaks down, 
we use the matched-asymptotic method introduced in ifTSl to 
find the DSW's slowly varying cnoidal solution. This naturally 
leads to Whitham's equations, which were originally obtained 
by the method of averaging [24] and later by a perturbative 
approach [34]. Then we use the method in |[T9l to determine 
the small, decaying, oscillating solution to the left of the DSW, 
which matches into the DSW. 

3.1. Shock front 

We asymptotically approximate u to the DSW's right by sum- 
ming the Neumann series for in (|7j. We do this by finding 
the long-time asymptotic-approximation of D(£; t) right of the 
DSW. Then we use the Neumann series formed from the iter- 
ates K (0> (x, y; t) = -HO + y; t) and K (n) (x, y; t) = -Q(x + y;t)- 
§™Q.{y + z;t)K { "- x \x,z;t)dz. 



Far to the DSW's right, the contribution to Q from the reflec- 
tion coefficient dominates and [u(x, t) + be 2 ] is exponentially 
small. 

Near the DSW's right, the contribution to Q. from the trans- 
mission coefficient dominates: the contribution from A — 
gives 



ef(f/r+4c 2 )/e 



Vi 



16 V^[6c - £/(2cf)] 3/2 



[H 2 (0)r 3/2 + O(r 512 )], 



where Hj(A t ) = [dJ\T(A,A r (A);0)\ 2 /dAj] A=A ,. 

The terms in the Neumann series become disordered when 
[x + 2c 2 1 + 3e/(4c) log(6c 2 f - x)] =0(1), which is at the DSW's 
right edge (cf. asymptotic principles discussed in 11351 ). The 
Neumann series can be summed, and we find that 



u(x, t) 6c + 12c sech 2 



where ft = £/(2c)log{327r 1/z /[H 2 (0)c l/z s i/z ]} and 
3e 

£ = -x- 2c 2 t log(6c 2 f - x) + Ai(x/t)r l + 

4c 



(8) 



(9) 



(We omit A i due to length.) This provides the boundary condi- 
tion on the DSW's right. 

This procedure gives the DSW's phase, ft. This phase only 
depends on //2(C)) (since Ho(Q) = Hi(0) = 0); the equivalent 
phase term lfl8l Eq. (2.25c)] for vanishing data is determined 
by r"(0) - [r'(0)] 2 /r(0), where r is the corresponding reflection 
coefficient. The equivalent phase term in the shock solution 
associated with Burgers' equation also depends on the initial 



data in a similar way (see 54.1 and |36|). 



3.2. DSW 

We find the DSW using matched asymptotics analogous to 
[ 18 1. First we make the variable change u(x, t) = -6c 2 + g((,, i) 
in ([TJ, based on (|S). Then we introduce the slow-variables Z = 
6£ and T = St (where 6 = 0(r x )) to get 



e 2 g((( + 88c ~ 4 ° 2 8( ~ 8t 



' 3e(3g 2 g^ + gg^ - 12c 2 g f ) ' 
" 4c(8c 2 r + Z) 



(10) 



To leading order, ( fTO) has the special solution (which can be 
found using the methods in [37]) 

g<£, t) ~ 4c 2 - V + 4eV(l - 2k 2 ) 

+ \2k 2 E 2 K 2 cn 2 - ft - Vi), k] , (11) 



where cn(z, k) is the Jacobian elliptic 'cosine' (see 113811 ). Here, 
k, k, and V are arbitrary constants when the right-hand-side of 
( [Tol l i s neglected but vary slowly in general. In the special case 
k = 1, k = c/s, and 7 = 0, g(£, t) = 12c 2 sech 2 [c(f - ft)/e] and 
exactly matches ([8]). 

As in 1 34], we use the multiple-scales method to determine 
k, k, and V, which vary with the slow-variables Z and T. This 
leads to three conservation laws, which determine k, k, and V. 



3 



First we introduce the rapid-variable 6((, t) with 6^ = k(Z, T) 
and 6, = -a>(Z, T) = -kV; this leads to the compatibility condi- 
tion (6^), = (6>,)^ or kj + a>z = — which is a conservation law. 
We then use d, = -wdg + 58t and d( = Kdg + Sdz to transform 
( fTO] ); then we expand g(0,Z,T) = g o (0,Z,T) + 6gi(0,Z,T) + 
irg 2 (6, Z,T) + ■ ■ ■ and group terms in like powers of 8. The 
solution of the 0(1) equation is 



Z, T) = a(Z, T) + b(Z, T) cn 2 {2K(6 - O ), k(Z, T)}, 



where K = K(k(Z, T)) is the complete elliptic integral of the 
first kind, k 2 = b/[48s 2 k 2 K 2 ], and a = 4c 2 - V-2b/3 + b/(3k 2 ). 
Then we enforce the periodicity of go(0, Z, T) in 9 to eliminate 
secular terms (that is, terms that grow arbitrarily large); this 
gives the other two conservation laws, which we omit due to 
length. These three conservation laws determine b, k, and V. 
If we make the variable change b/k 2 = 2(rj, - r{), k 2 = (r 2 - 



ri)/(r-) - n), and V = 4c 2 
gives the diagonal system 



(r\ + r 2 + r-})/3, then simplifying 



— + v i (r 1 ,r 2 ,r 3 ) — = 0, 
oT aZ 



i= 1,2,3, 



(12) 



where vi = V + bK/[3(K - E)], v 2 = V + b(l -k 2 )K/[3(E-(l - 
k 2 )K)}, v 3 = V- b(l - k 2 )K/(3k 2 E), and 



Z, T) = n - n + r 3 + 2(r 2 - n) cn 2 (2K(8 - 6 ), k) . 



Whitham first found this diagonal system in [24J (see also 

H3). 

For large-time, the solution tends to a self-similar solution: 
that is, r,- = r,-(y) with^ = Z/T = (It. The boundary conditions 
are satisfied when r\ — and = 6c 2 . So (jT2]» reduces to 

(v2 -x)Aix) = 0, 

which V2 = x satisfies. The numerical solution of this implicit 
equation for r 2 is plotted in Fig. [3] We can directly compute the 
right and left speed of the DSW: At the DSW's right, we take 

the limit r 2 — > and get that v 2 — > or x 2c 2 1 — the speed 

of the soliton train. At the DSW's left, we take the limit r 2 — > r\ 
and get that v 2 — > 10c 2 or x 12c 2 f. 



6 c 1 





(b) 
t/£=200 



Figure 3: The value of n(x) found numerically for < x < 10c 2 , where 
X = Clt. For comparison, we include —n(x) as a dashed line and a numerical 
simulation of u(x, t) in gray (inside the envelope of ri and — r^i) for a single-step 

at (a) t/e = 10 and (b) t/s = 200. Note that^ = corresponds to x 2c 2 ? 

and^ = 10c 2 to a: 12c 2 f. 



3.3. Trailing edge 

On the DSW's left, there is a decaying, slowly varying, 
oscillatory similarity-solution (as there is with vanishing data 



ED ED): 

x l/4 

u(x, t) ~ 2A cos(0) - 

where X = -x/(3t), t = 3t, 



A 2 (l - cos 26») 



+ 0(T- i/Z ), 



3 18 T T 



and A and 6> depend on the scattering data. 

We could, in principle, sum the Neumann series from the 
GLM integral equation formulated from -co to x in the long- 
time limit to find A and 6>o- Indeed, the GLM integral equation 
formulated from -oo to x has the same form for both step-like 
and vanishing data. 

Instead, we use the method in |19| to find A and 6\). Using 
this oscillatory similarity-solution in ( f3a] >, we asymptotically 
solve for the eigenfunction <p with boundary-values <p — > e^' A ^ E 
as x -> -oo and (f> -> ae~ iAxle - + be iA(x+&A2,)/e as x -> -12c 2 f. 
Here, a and b can be related to a and b through either the GLM 
integral equation (formulated from -oo to x) or the asymptotic 
forms of u for -12c 2 1 <K x <K oo. This leads to a matched- 
asymptotics problem with three regions — the resonance or 
turning-point region, where X ~ 4A 2 , is solved in terms of 
parabolic cylinder functions. The result is that 

9e 



A 2 (X)~--^log(l-|fl(Vx/2)| ) 



and 



7 = ? -arg{^,0)}-arg F 1-^— 



c 2 A 2 (4c 2 ) 
9sA 2 



log 



\c + A) 



6s 



log 2 



where A = y/X/Z, f(A,t) = b(A,t)la(A), and \f(A,t)\ = \R(A,t)\. 
The limits x -> -12c 2 f and r 2 ~ 2(10c 2 -^)/3 ~ 2AX 1/4 r-' /2 
yield u ~ 2 V2c/(3f) cos(16c 3 f/e), which matches the DSW at 
its left edge (see [36] for details). 

4. Discussion 

4.1. Comparison with Burgers' equation 

Burgers' equation (w, + ww x -vw xx — 0, v > 0) is the leading- 
order asymptotic equation for VSWs. If we take initial data that 
go rapidly to the boundary conditions lim^_<x, w(x, t) — and 
limj^ +00 w(x, t) = -h 2 , then the long-time asymptotic solution 

is w(x,t) (h 2 /2){l+tanh[h 2 (x-x + h 2 t/2)/(4v)]}, where x 

is a real constant that depends on the initial data [36|. So well- 
separated step data go to a single shock wave in the large-time 
limit for both Burgers' and the KdV equation. For both, the 
boundary data determine its form and the initial data determine 
its location. Unlike with Burgers' equation, the solution of the 
KdV equation with step-like data can also have a finite number 
of solitons, which move to the DSW's right in the long-time 
limit. 
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4.2. Comparison with the linear KdV equation 

The large-time asymptotic solution of ([l} differs significantly 
from the linear problem (5, + e 2 u xxx = 0). While the nonlinear 
problem has a central DS W region with strong nonlinearity over 
\x\ — O(t), the linear problem's middle region is only over \x\ = 
0(t 1 ^). Indeed, the solution to the linear problem in the middle 
m Ai(?7 ) drj , where Ai(x) is the Airy 
function, tj = x/(3s 2 t) 1 ^, and C/rj is the Fourier transform of 
U x (x, 0). 

4.3. Comparison with vanishing boundary conditions 

The large-time asymptotic solution of ([TJ for step-like data 
(c + 0) is also significantly different from that for vanishing 
data (c = 0). The collisionless shock (region III in Fig. [T]?), 
which is analogous to the DSW, has width <9[f I/3 (log f) 2/3 ]; the 
DSW (region B in Fig. [lj) has width O(t). From Ablowitz 
and Segur [18], the long-time asymptotic solution for vanish- 
ing data have four regions (Fig. [TJi) — all of which decay in 
time: an exponentially small solution in region I, x > O(t); a 
growing similarity-solution in region II, |jc| < 0(t 1 ^ 3 ), which 
is related to Painleve II; a collisionless shock in region III, 
(-x) = 0[t '^(log t) 2 ^ 3 ], which is a slowly varying cnoidal wave 
analogous to a DSW; and an oscillatory similarity-solution in 
region IV, (— x) > 0(f), which has the same form as the solution 
inf3] 

4.4. Conclusion 

In this letter, we show that general, step-like initial data 
tend to a single-phase DSW for large-time. Therefore, well- 
separated, multi-step initial data eventually form a single-phase 
DSW despite having multiphase dynamics at intermediate times 
(as indicated in Fig. [2]). The asymptotic solution of the KdV 
equation for general, step-like data is new. The details of our 
calculations will be given in a separate paper; we anticipate that 
they can be applied to other integrable nonlinear PDEs with 
general, step-like data. 
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